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ABSTRACT 

The lactose operon of Escherichia coii is a paradigm 
system for quantitative understanding of gene regu- 
lation in prokaryotes. Yet, none of the many math- 
ematical models built so far to study the dynamics 
of this system considered the fact that the Lac re- 
pressor regulates its own transcription by forming 
a transcriptional roadblock at the 03 operator 
site. Here we study the effect of autoregulation on 
intracellular Lac! levels and also show that cAMP- 
CRP binding does not affect the efficiency of 
autoregulation. We built a mathematical model to 
study the role of Lad autoregulation in the lactose 
utilization system. Previously, it has been argued 
that negative autoregulation can significantly 
reduce noise as well as increase the speed of 
response. We show that the particular molecular 
mechanism, a transcriptional roadblock, used to 
achieve self-repression in the lac system does 
neither. Instead, Lac! autoregulation balances two 
opposing states, one that allows quicker response 
to smaller pulses of external lactose, and the other 
that minimizes production costs in the absence of 
lactose. 

INTRODUCTION 

Bacteria sense a wide array of signals (minerals, nutrients, 
stress signals, etc.). A large class of cellular response 
systems regulates the flux and concentration of small mol- 
ecules by controlling transport and metabolism pathways 



via two feedback loops connected by a common transcrip- 
tion regulatory protein that senses the intracellular con- 
centration of the small molecule (1,2). In fact, almost half 
of the transcriptional regulators in Escherichia coii are 
directly regulated by a small molecule (3). The prototypic 
example of such a control system is the lac operon, which 
has been a paradigm of gene regulation. In E. coii, the lac 
operon contains genes encoding the lactose transporter 
(LacY) and the enzyme for lactose degradation (LacZ), 
therefore the lactose repressor (Lad) regulates the trans- 
port and metabohsni pathways simultaneously (4). The 
lad gene is present just upstream of the lac operon, and 
in fact there are three operator sites where the Lad 
tetramer can bind and affect transcription (5). The Lad 
tetramer contains two identical dimers, connected at their 
C-terminal region. Each dinier in the tetrameric structure 
has an N-terminal helix-turn-helix DNA-binding domain 
(6). The structure of the lac system is shown schematically 
in Figure 1. The main operator is 01, the strongest of the 
three operator sites. Lad binding to 01 represses tran- 
scription of the lac operon but leaves the expression of 
the lad gene unchanged. The binding of Lad to 01 in- 
creases its probabihty to bind via DNA looping to 02 or 
OJ, which are weaker operators (7,8). When bound to 01 
and (92, transcription of the lac operon is repressed, while 
Lad continues to be produced. However, when 01 and 
03 are bound, not only is the lac operon repressed, but the 
production of Lad is also prevented (9). In this state, 
transcription of lad occurs but only a truncated transcript 
is produced, which is in turn subject to SsrA-mediated 
tagging and subsequent proteolysis of the truncated 
protein produced (9). While there is experimental 
evidence for Lad autoregulation (9,10), this feature of 
the network is ignored by the available mathematical 
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Figure 1. Schematic structure of the lac genes and control regions on the E.coli chromosome and in the pSEMJl plasmid, used as a template for 
in vilro transcription. Arrows represent promoters, smaller dark gray and light gray boxes represent Lad and cAMP-CRP binding sites, respectively. 
The lac region in the pSEMJl plasmid is flanked by transcription terminators (T). The Piad - 03 distance in pSEMJl is 1022 bp shorter than on the 
chromosome but the 03-02 region has the same sequence as the natural chromosomal gene. Transcripts initiated at P/„^i can be roadblocked by Lad 
at 03 (RBI) or at Ol (RB2). The binding states considered in the model are shown in the bottom. Relative probabilities of the Lad-bound states 
depend on the active Lad concentration (/*) and on the associated binding energies (s, = fiAC/tr ^ 0.6nM"', x 13.8nM"', Ej 28.9nM"'. 
The figure is not drawn to scale. 



models (11-15). Previous studies suggested that negative 
autoregulation in regulatory networks can significantly 
reduce noise (16) and increase the speed of response 
(17). In this work, we study the effect of autoregulation 
on intracellular Lad concentration and build a stochastic 
model of the lactose utilization system to explore the role 
of Lad autoregulation. We compare the natural lac 
system with two hypothetical controls, where Lad is 
produced at a constant low or at a constant high level, 
which correspond to the estimated autoregulated and fully 
expressed Lad levels, respectively. We show that the 
mechanism of Lad autoregulation neither reduces noise 
nor increases the speed of response. However, we find that 
the autoregulated system has a larger dynamic range and 
performs more economically than the constitutive systems. 

MATERIALS AND METHODS 

Plasmid construction 

The pSEMJl plasmid, used as a template for in vitro tran- 
scription, was created by inserting the _p£^^ promoter 
region from plasmid pTYBl (NEB) and the Piacz 
promoter region {03-02, nt 365820 365101) froin 
E. coli MG1655 (GeneBank: NC_000913.2) between the 
EcoKl and Pst\ sites of plasmid pSEM2008. The Pf^^^.^ 
promoter region was polymerase chain reaction (PGR) 
amplified using the primers ATATATGAATTCGAATG 
TTGACAAACCTTTCGCGGTATGGCATGATAGC 
and ATATATCTCGAGATTCACCACGCTGAATTGA 



CTCTCTTC to replace the original -35 promoter 
element with the TTGACA consensus sequence (the re- 
sulting enhanced promoter is termed i^^^./). The Pjacz 
promoter region was amplified using the primers ATAT 
ATCTCGAGCAACTCTCTCAGGGCCAGGCGGTG 
AAGGGC and ATATATCTGCAG AATAATTCGCGT 
CTGGCCTTCCTGTAGCCAGC. The Pj;^,^ PGR 
fragment was cut with EcoKl and Xho\, and the Piacz 
PGR fragment was cut with Pst\ and Xhol. The two frag- 
ments were inserted between the EcoRl and Pst\ sites of 
plasinid pSEM2008 by a three-piece hgation. The 
pSEM2008 plasmid was obtained by inserting the DNA 
fragment containing the rrnBTlTl terminators (nt 4559 
4141) froin pKK223-3 (Pharmacia, GeneBank 
M77749) between the Kpnl and EcoKl sites of 
pSA850 (18). 

To create plasmid pSEM1068 for the expression of the 
Hise-tagged dimeric Lad protein (lacking the last 16 
amino acids), the lad gene was PGR amplified using the 
primers AAAAGGTAGGAAAAGGTTTGGGGGTATG 
GGTGAT and AAAAGAATTGAAGGGAA GGAGGT 
GGATGGGGGAAG, the amplified DNA fragment was 
digested with Nhel and iscoRI, and inserted into the 
pSEM1026 vector (19) between the Nhel and EcoRR 
sites. The sequence of the ampUfied region was verified. 

Protein purification 

The Hisg-tagged dimeric Lad protein was expressed in 
E. coli strain Top 10 bearing pSEM1068 and purified 
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using the protocol described previously for the purifica- 
tion of Hisg-tagged GalR (19). CRP was purified as 
described by Ryu et al. (20). The wild-type (WT) Lacl 
protein was a kind gift from Maxim Sukhodolets. 

In vitro transcription and quantitation 

Reactions were performed on supercoiled pSEMJl 
plasmid DNA as described previously (21). Lad and 
CRP were used at the concentrations shown in Figure 2, 
when present. The RNA bands were quantified using the 
Storm 860 Phosphorlmager (GE Healthcare). The lengths 
of the roadblocked transcripts were estimated based on 
the migration distances of transcripts with known lengths. 

Western blot and quantitation 

Protein samples were loaded on a 10% Bis-Tris gel as 
follows: 10'°, 5 X 10'°, lO", 2.5 x lO" Lacl repressor mol- 
ecules mixed with cell extracts obtained from 2x10^ 
E. coli MC4100 {MacJ) cells, extracts of 2 x 10^ 
MG1655 cells that were grown in the presence and 
absence of isopropyl P-D-thiogalactopyranoside (IPTG; 
1 mM), respectively, and extracts obtained from 2x 10** 
E. coli MC4100 cells. 

After separation, the proteins were transferred to an 
Immobilon-P PVDF Membrane of 0.45 |im pore size 
(Millipore), then were blocked overnight with 5% nonfat 
dry milk in PBST (50 mM phosphate buffer at pH 7.4, 
650 mM NaCl and 0.1% Tween 20). The blot was then 
incubated with primary antibody overnight at 4°C. It was 
then incubated with peroxidase-conjugated antibody and 
developed with Supersignal West Pico kit (Thermo 
Scientific). The dilutions of the antibodies were anti-LacI 
antibody (1:1000, Millipore), peroxidase-conjugated anti- 
mouse antibody (1:2000, Sigma A2554). Band intensities 
were quantified and background corrected. The bottom 
signal, which is also present in the E. coli MC4100 
extract, was used as an internal control. 



Lacl (nM) 0 18 120 12 0 18 120 12 
CRP(nM)0 0 0 0 50 50 50 50 




Mathematical model 

The dynamical variables we keep track of in our model are 
the concentrations of internal lactose (L), internal 
allolactose {A), Lad mRNA (/,„), Lacl tetramers (/), 
LacY permeases (Y) and the LacZ enzymes (Z). The de- 
terministic differential equations that model the dynamics 
of these variables are described below: 



dA ^ L ^ A 

— v.Z v-Z 

dt - L+K - A+K 



Iv-Z- 



L+K 



(1) 
(2) 



The first term on the right side of equation (1) represents 
the import of external lactose (Lpxr) by LacY, for which we 
have chosen a Michaelis-Menten form where v,. is the 
maximum rate of import per LacY molecule, and the con- 
stant Kgy-t is the Lpvr concentration at which the import 
rate per permease is half of its maximum value. The 
second term is similar and represents the export of 
internal lactose by LacY. The export has a different 
Michaehs constant, which is larger {XK^^,). The final 



RNA1 



108 
106 



1 



5 6 



8 



Figure 2. Effect of Lad and cAMP-CRP on transcription from the 
P'Lci promoter in vitro. Concentrations of WT Lacl (lanes 2 and 6), 
dimeric Lacl (lanes 3, 4, 7 and 8) and CRP (lanes 5-8) are shown on 
top. cAMP was present at 100|.iM concentration in all reactions. Two 
RNA species appeared in the presence of WT Lacl as a result of road- 
block termination (labeled RBI and RB2). The RNAl transcripts (106 
and 108 nt) were used as internal controls between lanes. The lengths of 
RBI and RB2 were estimated based on the migration distances of the 
RNAl transcripts and of the transcripts initiated at promoters Piacz 
(506 nt), Pb (322 nt) and Pc (404 nt) (36). 



term models both hydrolysis of internal lactose as well 
as its conversion to allolactose, catalyzed by LacZ. 
Again a Michaelis-Menten form was chosen for these re- 
actions, with each having a maximum rate of v- per LacZ. 
We further assume that both reactions have the same 
Michaehs constant K (22). The terms on the right side 
of equation (2) similarly model the production of 
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allolactose and its hydrolysis by LacZ with the same 
Michaehs constant K (23). We assume that the binding 
and unbinding of allolactose to Lad is fast, so that we 
can take the concentration of active Lad (i.e. unbound to 
aholactose) to be 



r 



1- 



■(*)' 



(3) 



where is the Michaehs constant of allolactose-LacI 
binding and h is the associated Hill coefficient. Active 
Lad can bind to three operators 01, 02 and 03. Lad 
tetramers bound to these operators can interact by the 
formation of DNA loops. If 01 is bound, with or 
without a DNA loop, we assume that transcription of 
the lac operon is completely blocked. If 03 is bound, we 
assume that transcription of lad is roadblocked. In prin- 
ciple, there are then 14 possible states each with their par- 
ticular combination of transcriptional repression of the 
lad gene and the lac operon. In our model, we wiU only 
allow four states: (i) all operators are unbound, (h) 01 is 
bound, (iii) a DNA loop is formed between 01 and 02 
and (iv) a DNA loop is formed between 01 and 03. This 
is because the other bound states have significantly lower 
energies, and hence lower probability of occurrence, than 
the three bound states we allow (7). Promoter activity 
levels can be given as the sum of the products of 
promoter activities and probabihties for all possible 
states (24). Based on Figure 1, we can write the activities 
of the lac operon and the lad gene as a function of the 
active Lad concentration (I*): 



lac activity {r) 



lad activityir) 



1 



l+(Si+£2+£3)/* 

l+(£l+S2)/* 
" l+(Si+£2+£3)/* 



(4) 



(5) 



The £s are related to the binding energies of the Lad- 
operator complexes, £,• — e'^'^il^T , in addition, the 
activity of the lac operon is controlled by the cyclic- 
AMP-CRP level, C: 



lac activity {T ,C) ■ 



1 a+£cC 

l+(£l+£2+£3)/* ^ l+£cC 



(6) 



The effect of cAMP-CRP is taken to be completely inde- 
pendent of the effect of Lad on the total activity, as we 
show in Figure 2. 

Using equations (5) and (6) we can write the differential 
equations for the relevant mRNA and protein levels: 



dim 

dt 



l+(El+£2)/* 



_1+(Ei+£2+£3)/*. 



YmJm 



di , 

dY , 



dt 



a+EcC 



[l+(£l+£2+£3)/*][l+Scq 



yY 



(7) 

(8) 
(9) 



dZ 
dt 



a+SrC 



■[l+(£l+£2+£3)/*][l+£cC| 



yZ 



(10) 



The /cs are parameters that set the maximal rates of tran- 
scription and translation, and the ys set the dilution and 
degradation rates of the proteins and lad mRNA. We 
have chosen to model the lad niRNA explicitly but not 
the lac operon mRNA because the former is produced at a 
sufficiently low rate (1 transcript per generation) (25,26) to 
produce significant fluctuations in the Lad levels, whereas 
the lac operon mRNA is produced at a high enough rate 
to have little effect on the fluctuations of LacY and LacZ 
levels. 

For these deterministic differential equations, there is 
no need to have a separate equation for dZ/dt. Instead 
Z can simply be calculated from Y because 
Z(0 = Z(0)+|;[F(0 - F(0)]. We simplify the modd one 
step further by assuming that the transport of lactose 
and hydrolysis of lactose/allolactose take place much 
faster than the processes of transcription and translation. 
Then we can assume that L and A concentrations are 
always in quasi-equihbrium. Setting dL/dt = dA/dt = 0 
makes the lactose and aholactose concentrations both 
equal to the physically sensible (i.e. real and non- 
negative) solution of the following: 



C2L +ciL — Co = 0 
where 

Co = XKexlKVy YLexl, 

c\ = VyYK{l - L„,) - VyYXKextLex&'^v^ZXKi.xt, 



Ci 



VyYil - L^J+2v,Z, 



Lext = Lext / {Kcxt~^ Lext)- 

So, 



1 

2^2 



-cr 



Cf+4C2C0 



(11) 



(12) 
(13) 
(14) 
(15) 



(16) 



and the rest of the equations remain unchanged. Thus, 
equations 3, 7-10, and 12-16 are used for the deterministic 
simulations. 

In the stochastic simulations, we keep track of the 
actual number of molecules of Lad tetramers. Lad 
mRNA, LacY and LacZ tetramers. The net production 
and degradation rates, from the equations of the deter- 
ministic model, expressed in appropriate units, can be 
treated as probabihties per unit time for the production 
and degradation of each species. We use the Gillespie al- 
gorithm (27) to determine, from these probabihties per 
unit time, the time at which the next production or deg- 
radation will happen and which species it will affect. We 
then accordingly increase or decrease the number of that 
species, recalculate the probabilities per unit time of pro- 
duction and degradation and repeat. This gives us a time 
series of the number of Lad tetramers, Lad mRNA, 
LacY and LacZ tetramer molecules as a function of 
time. Figures 5 and 6, and all statistics from them, were 
obtained from an ensemble of such stochastic simulations. 
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The number of parameters used in the model could be 
reduced by non-dimensionalizing the equations but this is 
not necessary here because we can fix most of the param- 
eter values. Only the external lactose concentration 
is varied in the simulations, the remaining parameters are 
always kept fixed. 

Parameter values: 

The value of y = 0.0087 min"' = ln(2)/(80min) was 
chosen assuming a doubling time of 80 min (28) and no 
specific degradation of the proteins. The values used for 
the half-saturation constant for active transport by LacY 
(Ar„/) and the active transport turnover number (v,,) were 
0.27 mM and 48/s, respectively (29). 

The £s (s/ = e^*^'/''^) can be determined from the fold- 
repression of promoter activities in different conditions. In 
the WT cell with autoregulated level of Lad tetramers 
(30 nM), the repression of the lac operon activity is 
1300-fold, and in a cell that has only the Ol operator, 
the repression is 18-fold (7). Assuming that binding to 
the operator sites is sufficiently strong (to be precise, 
assuming /*(£i+£2) >> 1), this sets « 18/30nM~', and 
(£i + e2 + S3) Si 1300/30nM~'. Further, the experimental 
results presented in this article set the ratio between the 
autoregulated and the fully expressed levels of Lad to be 
1/3. Therefore £i+S2 = £3/2, i.e. £2 = 415.33/30 nM"' and 
£3 = 866.66/30 nM"'. This choice results in 1300-fold re- 
pression in the deterministic simulations but a lower mean 
repression is observed in the stochastic simulations 
because of the noise in Lad levels. The repression level 
can be increased by increasing the Lad-operator binding 
energies or by increasing the number of lad mRNAs 
produced per cell generation. We have tested these 
possibilities and the conclusions reported here were not 
affected. 

We study our model only in conditions where the 
cAMP-CRP level is fixed and close to saturation, so we 
take the factor (q' + £cC)/(1 + e^C) = 0.9. 

For the allolactose-Lacl interaction, we use a Hill co- 
efficient of 2 and = 1 ^M (12,30), while for LacZ- 
mediated hydrolysis of lactose and allolactose, we use 
A'=1.4|rM (22). The maximal transcription rate of 
the lad mRNA (/c,.) is set to l/80nM/min because the 
lad gene is transcribed approximately once per cell 
generation on average (25,26). The half-hfe of the lad 
mRNA is ~3.8min (31), therefore y,„ = In 2/ 
3.8 min = 0.1824/min. The maximal transcription rates 
of lacY (k^) and lacZ (k-) were chosen to be 90 and 
lOOnM/min, respectively, to obtain ~10(iM LacY and 
slightly higher LacZ tetramer concentration when Lad 
is inactivated (32,33). 

The combined conversion rate of lactose by LacZ (hy- 
drolysis plus conversion to allolactose) is 3600/min (11), 
therefore vz = 1800 min"'. 

The maximal translation rate of the lad mRNA (ki) was 
chosen to be 90 jjm/kc per min to obtain 90 nM Lad 
tetramers/cell at saturating intracellular allolactose 
concentration. 

This leaves one remaining parameter, 1, which sets the 
rate of LacY-mediated lactose export. To fix this, we need 



one additional constraint. We use the approximation that 
~3 bilhon glucose molecules are needed to generate a new 
cell (34). We found that with k = 750 (dimensionless), 
with 5mM external lactose concentration the cell, in our 
model, metabolizes x 2 bilhon lactose molecules in 
80 min, which must provide sufficient resources to 
generate a new cell. 



RESULTS 

In vitro transcription pattern of the lad mRNA 

In vivo studies of the distribution of different lad mRNA 
species in rifampicin-treated cells showed that in the 
presence of Lad, transcription is blocked at the 03 and 
01 operators with similar efficiencies (10). We constructed 
a plasmid DNA template that allows us to study how 
termination of transcription is controlled in vitro. 
The plasmid constructed (pSEMJl) contained an 
enhanced P,,,,., promoter and the WT lac control region 
including both the 03 and 02 (Figure 1). The reason for 
using a control region spanning the whole 03-02 region 
instead of the 190-bp control region containing only 03 
and 01 (10) was that a Lad tetramer can bind either 01 
and 03 or 01 and 02, and only the 03-01 loop can 
regulate lad transcription. So far it is unclear how the 
choice between 02 and 03 binding is made by an 07- 
bound Lad tetramer. Lad bound to the lac control 
region efficiently separates the lad and lacZYA transcrip- 
tion units; >90% of lad transcription is blocked when 
both auxiliary operators are present (35). We performed 
in vitro transcription assays using supercoiled pSEMJl 
plasmid DNA to study the effect of Lad tetramers. Lad 
dimers and cAMP-CRP on the elongation of the lad 
mRNA (Figure 2). In the presence of WT Lad, two 
RNA species appeared as a result of premature termin- 
ation (labeled RBI and RB2). However, in the presence of 
Lad dimers only a faint band appeared at the position of 
the shorter transcript (RBI), while the longer transcript 
(RB2) had similar intensity to the one obtained with WT 
Lad, suggesting that RBI corresponds to a transcript 
terminated by Lad occupying 03, and RB2 is terminated 
owing to Lad binding to 01. Length estimation of the 
roadblocked transcripts suggest that RBI and RB2 cor- 
respond to the in vivo observed mRNAs with endpoints 
Iir and II (10), respectively. The effect of cAMP-CRP on 
03-01 loop formation is ambiguous in the hterature 
(37,38). In our assays, the presence of both WT Lad 
and cAMP-CRP resulted in the same pattern obtained 
with WT Lad only (Figure 2, lane 2 versus lane 6). The 
amount of RB2 transcript obtained in the presence of Lad 
dimers was shghtly reduced in the presence of cAMP-CRP 
(Figure 2, lanes 3 and 4 versus lanes 7 and 8). Therefore 
we concluded that even if cAMP-CRP is bound to the 03- 
01 loop, its presence does not affect the quantity of the 
RBI and RB2 roadblock products, and it does not act as 
an additional roadblock. Therefore, in the mathematical 
model, the effect of cAMP-CRP is taken to be completely 
independent of the effect of Lad. 
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Figure 3. The effect of autoregulation on in vivo Lad levels. 
Escherichia coli MG1655 cells were grown in the absence (lanes 5 
and 7) and in the presence (lanes 6 and 8) of 1 mM IPTG and 
proteins of 2 x 10** cells were separated by sodium dodecyl sulpha te- 
polyacrylamide gel electrophoresis and analyzed by western blotting, 
using an anti-LacI antibody. Purified Lad repressor molecules were 
mixed with proteins from 2x 10** cells of E.coli MC4100 (A/of/) (41) 
and loaded as controls (lanes 1^, 2.5 x lO", lO", 5x 10'°, and 10'° 
Lad monomers, respectively. Proteins from 2 X lO'* cells of E.coli 
MC4100 were loaded in lane 9. The amount of Lad in cells was 
estimated by quantifying the signals obtained. 



The effect of autoregulation on in vivo Lad levels 

Most of the current models assume that cells contain ~10 
Lacl tetramers, based on the estimate of Gilbert and 
Miiller-Hill (39). A more recent work predicted 8.8 tetra- 
mers on average in HG104 {AlacZYA) cells as a lower 
bound on the actual number of proteins in vivo (40). We 
measured the average Lad content of cells in the presence 
and absence of IPTG, corresponding to fully expressed 
and autoregulated levels respectively, using western 
blotting of known amount of cells (Figure 3). Based on 
the experiment, we estimate that ~ 150-1 80 monomers are 
present on average per cell in the presence of IPTG, while 
only 55-65 monomers/cell can be found in the absence of 
IPTG. Based on our estimate, the upper bound on the 
number of Lad tetramers is ~40 molecules on average 
in the fully induced cells and 15 molecules in the absence 
of inducer. The latter falls between the estimate of 8.8 
tetramers in HG104 cells and the theoretical calculation 
of 30 nM Lad tetramers by Santillan and Mackey (11). 

Development of the mathematical model for the lac system 

We have developed a mathematical model to study the 
effects of Lad autoregulation. The model is described in 
detail in the 'Materials and Methods' section. Using this 
model, we compared the steady state and dynamic behav- 
iors of the natural lac system, where Lad expression is 
autoregulated, with two hypothetical controls where 
Lad is produced at a constant low or at a constant high 
level, which correspond to the estimated autoregulated 
and fully expressed Lad levels, respectively. Because of 
the uncertainty of intracellular concentration of Lad 
tetramers, we have performed two sets of each computa- 
tion, reflecting the higher and lower estimates found in the 
hterature. In the first set, the average autoregulated Lad 
level was 30 nM and the fully expressed level was 90 nM, 
while in the second set these values were 10 and 30 nM, 
respectively. 

Steady-state simulations 

In case of the lac system the input dynamic range can be 
defined as the extracellular lactose concentration (input) 
interval over which the average level of LacY transporter 
(output) changes significantly. Previously, systems 
controlled by negatively autoregulated regulators were 
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Figure 4. Average number of LacY molecules per cell as a function of 
extracellular lactose concentration, obtained by numerically solving the 
deterministic differential equations that define the model (see 'Materials 
and Methods' section). The black curve represents the WT system, 
while the dark gray and light gray curves represent the systems with 
non-regulated Lad expression, corresponding to 30 nM (fixed-low) and 
90 nM (fixed-high) Lad levels, respectively. For large external lactose 
concentrations, the LacY number approaches ~9400 molecules (we 
assume 1 nM is approximately one molecule per cell). Dynamic range 
(external lactose level at which 90% of this level is reached divided by 
external lactose level at which 10% of this level is reached) is ~14.4 for 
WT, 9.2 for fixed-low and 9.4 for fixed-high systems. 



found to have a larger input dynamic range and a more 
Hnear dose-response compared to similar systems 
regulated by constitutively expressed regulators (42,43). 
To explore the effect of Lad autoregulation on the 
input dynamic range of the system we have computed 
the average level of the LacY lactose transporter per 
cell, at different extracellular lactose levels, using deter- 
ministic simulations (Figure 4). In defining the input 
dynamic range as the ratio of extracellular lactose levels 
at which the system shows 90 and 10% of its maximal 
output, we follow Goldbeter and Koshland (44). The 
plot shows that the WT system has ~50% larger 
dynamic range compared with the systems where Lacl is 
present at constant low (30 nM) or high (90 nM) levels. A 
quantitatively similar increase in the dynamic range was 
obtained when Lacl concentrations ranged from 10 to 
30 nM instead. 

The lacl gene is transcribed from a weak promoter re- 
sulting in about one new lacl mRNA per cell generation 
(25,26). Due to the noise originating in stochastic intracel- 
lular processes and the low number of short-lived lacl 
mRNAs, repressor levels fluctuate with time within each 
cell and differ among isogenic cells (45). Negative 
autoregulatory feedback loops in gene circuits have been 
shown to limit the range over which the concentrations of 
network components fluctuate (16). We performed sto- 
chastic simulations at zero extracellular lactose concentra- 
tion to compare Lad and LacY distributions in the 
constitutive and autoregulated systems (Figure 5). 
Interestingly, similar Lad levels were found in the WT 
system (mean ± SD = 33.1 ± 43.7) and in the system 
having constant low level of Lad (30 ± 42.7). Thus, 
autoregulation does not seem to reduce noise in Lad 
levels compared with the constitutive low system. 



Nucleic Acids Research, 2013, Vol. 41, No. 13 6387 



Nevertheless, the LacY level in the model is much higher 
for the system where LacI is expressed at a constant low 
level (1306 ± 2566) compared with the WT system 
(491.8 ± 1271). This is because probabiUty of having 
zero intracellular LacI was about three times lower in 
the WT system (Figure 5 inset), which influences the 
average LacY levels significantly. We also computed the 
protein levels for the system where LacI is expressed at a 
constant high level (90 nM). Because of its higher level, the 
intracellular LacI level fluctuates less (90 ± 74 nM), and 




0 2000 4000 6000 8000 10000 

Number of LacY molecules 




0 100 200 300 400 500 600 



Number of LacI tetramers 

Figure 5. Distribution of the numbers of LacY (top) and LacI 
(bottom) molecules from stochastic simulations with zero external 
lactose (duration of 1 1 million min). The first 1 million min were dis- 
carded before making these distributions to eliminate any transients. 
The black curve represents the WT system, while the dark gray and 
light gray curves represent the fixed-low and fixed-high systems, re- 
spectively. Table 1 lists the mean and standard deviations for each 
distribution. 



results in a more successful repression of the lac operon. 
As a consequence of stronger repression, the LacY level 
was lower and the variability of the LacY level was found 
to be substantially higher (29.9 ± 301.6 nM). Simflar 
results were obtained in a second set of simulations 
where the LacI levels were in the range of 10-30 nM 
(Table 1). 

Dynamic simulations 

Escherichia coli cells need to optimize their gene expres- 
sion pattern in environments where the quality and 
amount of carbon sources fluctuate, most likely in an un- 
predictable fashion (46). To test how fast the WT system 
responds to changes in extracellular lactose levels 
coinpared with the systeins having constitutive LacI ex- 
pression, we performed stochastic simulations where the 
external lactose level was changed from 0 to 5mM, and 
later back to OnM. We recorded the turn-on and turn-off 
times, which are defined as the time taken to reach 95 and 
5% of the maximal LacY levels, respectively (Figure 6). 
We find that the system expressing LacI constitutively at a 
low level has a longer turn-off time on average, and higher 
population heterogeneity in both turning on and off 
(Table 2). Furthermore, although the average turn-on 
times are similar, the system with fixed high LacI almost 
always takes much more than a ceU generation to turn on, 
whereas some ceUs in the WT and fixed-low systems turn 
on at times even much less than a cell generation. More 
precisely, at 299 min after the concentration of external 
lactose jumped from 0 to 5 mM, we found that 0 out of 
1000 cells with fixed high LacI had turned on, whereas 52 
and 139 cells out of 1000 had turned on in the WT and 
fixed low systems, respectively. Again, similar results were 
obtained in the second set of coinputations where LacI 
ranges from 10 to 30 nM. 

DISCUSSION 

Autoregulation is a common feature of sugar-specific tran- 
scription regulatory proteins in E. coli (47). Negative 
autoregulation typically reduces the rate of transcription 
initiation (48-50), however, in the lactose system of E. coli 
transcription elongation is inhibited (9). The lac system is 
intrinsically noisy because of the low probability of lad 
transcription (26) and because of the topology of the regu- 
latory elements, i.e. simultaneous transcriptional regula- 
tion of both the lad gene and the /acZYA operon by a 



Table 1. Mean ± standard deviations of the numbers of LacI and LacY molecules from stochastic simulations, with zero external lactose 
(Figure 5) 





LacI ran{ 


je: 30-90 molecules 




LacI 


range: 10-30 molecules 






WT (autoregulated) 


Constant low 


Constant high 


WT (autoregulated) 


Constant low 


Constant high 


LacI molecules 


33.1 ± 43.7 


30.0 ± 42.7 


90.0 ± 74.0 


11.9 ± 15.1 


10.1 ± 14.6 


30.1 ± 25.2 


LacY molecules 


491.8 ± 1271 


1306 ± 2566 


29.9 ± 301.6 


851.7 ± 1658 


2119 ± 3065 


113.4 ± 651.2 



The last three columns are for simulations where the rate of translation of LacI mRNA is one-third of the default value described in 'Materials and 
Methods' section, thereby LacI levels range from 10 to 30 tetramers (the binding energies were appropriately modified to have the same repression 
levels). 
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Figure 6. Distributions of turn-on time (left) and turn-off time (right) obtained from 1000 stochastic simulations where the external lactose con- 
centration was zero for t = Omin to t = 10000 min, then external lactose was fixed at 5mM for t = 10000 min to t = 20000min, and finally external 
lactose was again set to zero from t = 20 000 min to t = 30 000 min. The long times between changes of external lactose were chosen simply to allow 
enough time for the system to reach steady state before each change. Turn-on time for each simulation was the time after t = 10 000 min required for 
the LacY level to first reach 9025 molecules. The black curve represents the WT system, while the dark gray and light gray curves represent the fixed- 
low and fixed-high systems, respectively. Turn-off time for each simulation was the time after t = 20 000 min required for the LacY level to first reach 
475 molecules. Table 2 hsts the mean and standard deviation for each distribution. 



Table 2. Mean ± standard deviations of turn-on and turn-off times obtained in the simulations shown in Figure 6 





Lad ran{ 


It: 30-90 molecules 


Lad 


range: 10-30 molecules 




WT (autoregulated) 


Constant low 


Constant high 


WT (autoregulated) 


Constant low 


Constant high 


Turn-on time (minutes) 
Turn-off time (minutes) 


354.2 ± 41.7 
363.7 ± 74.4 


330.5 ± 85.4 
468.9 ± 246.0 


368.8 ± 36.0 
348.4 ± 33.0 


349.3 ± 44.2 
394.0 ± 128.6 


308.7 ± 103.7 
602.0 ± 392.5 


364.8 ± 38.3 
354.4 ± 50.5 



The last three columns are for simulations where the rate of translation of Lad mRNA is one-third of the default value described in 'Materials and 
Methods' section, thereby Lad levels range from 10 to 30 tetramers (the binding energies were appropriately modified to have the same repression 
levels). 



single Lad tetramer bound to 01 and 03. The intrinsic 
noise can generate heterogeneous expression of LacY and 
LacZ in the cells of E. coli populations in the absence of 
lactose (45). The lac system was also found to exhibit 
bistable behavior in the presence of low levels of non- 
metabolized inducers such as thiomethyl P-D-galactoside 
and IPTG but no bistabiHty was observed when lactose 
was used as an inducer (14,51). 

The capabihty of negative transcriptional 
autoregulation to reduce gene expression noise has been 
demonstrated both experimentally and theoretically 
(16,52). However, strong negative autoregulation was 
also reported to have the opposite effect, an increase in 
the protein variabihty sacrificed for reduced mRNA usage 
(53). Cell to cell heterogeneity can be either beneficial or 
disadvantageous, and therefore regulatory systems may 
evolve either to reduce or maintain it. We have addressed 
the effect of negative autoregulation of Lad on protein 
level variability in the absence of inducer theoretically, by 
building a mathematical model and performing stochastic 
simulations. Results of simulations suggest that the noise 
in Lad expression in the absence of inducer, quantified as 
the squared coefficient of variation (54), is not smaller in 
the WT (negatively autoregulated) system than in the con- 
stitutive systems expressing constant low or high levels of 
LacL Interestingly, even though the mean Lad level and 



noise is similar in the negatively autoregulated and consti- 
tutive low systems, the mean LacY level was found to be 
about two times higher in the constitutive system. This is 
because it is more probable to have zero Lad in a cell in 
the constitutive system. Cells having zero Lad express 
genes of the lac operon at a high level but can use 
lactose immediately when it becomes available. Negative 
autoregulation of Lad decreases the probability of having 
zero Lad in the cell, from which we can speculate that in 
natural habitats the cost associated with higher LacY ex- 
pression and slower turn off of the constitutive low system 
is higher than the potential benefit of fast lactose utiliza- 
tion in a fraction of the population. Unlike the constitu- 
tive low system, the constitutive high system is at least as 
economical as the WT system in the absence of lactose — 
although it expresses Lad at a higher level, expression of 
the lac operon is much lower in this system. However, cells 
with constitutive high Lad levels would perform worse 
than constitutive low or WT cells in an environment 
where these systems compete for small pools of lactose 
that appear rarely and intermittently. This is because a 
fraction of the populations having the constitutive low 
or WT systems can start to utilize the lactose source 
much earlier than the average turn-on time, while the con- 
stitutive high system lacks this opportunity. 
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Because the Lad level can be changed in the WT 
system, this system balances the two opposing states, 
one that allows quicker response to smaller pulses of 
external lactose, and the other that minimizes production 
costs in the absence of lactose. The resulting increased 
dynamic response range therefore enhances the overall 
performance of the autoregulated lac system. 
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